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The structural properties of single component fluids of hard hyperspheres in odd space dimen- 
sionalities d are studied with an analytical approximation method that generalizes the Rational 
Function Approximation earlier introduced in the study of hard-sphere fluids [S. B. Yuste and A. 
Santos, Phys. Rev. A 43, 5418 (1991)]. The theory makes use of the exact form of the radial dis- 
tribution function to first order in density and extends it to finite density by assuming a rational 
form for a function defined in Laplace space, the coefficients being determined by simple physical 
requirements. Fourier transform in terms of reverse Bessel polynomials constitute the mathematical 
framework of this approximation, from which an analytical expression for the static structure factor 
is obtained. In its most elementary form, the method recovers the solution of the Percus-Yevick 
closure to the Ornstein-Zernike equation for hyperspheres at odd dimension. The present formalism 
aUows one to go beyond by yielding solutions with thermodynamic consistency between the virial 
and compressibility routes to any desired equation of state. Excellent agreement with available 
computer simulation data at d = 5 and d = 7 is obtained. 

PACS numbers: 61.20.Gy, 61.20.Ne, 05.20.Jj, 51.30.+i 



I. INTRODUCTION 

Fluid systems made of hard bodies constitute simple 
models in which impenetrable particles interact solely 
through hard-core repulsions. Despite their simplicity, 
hard-sphere systems are commonly used as reference 
models to obtain accurate descriptions of real substances 
(simple fluids, colloidal suspensions, granular media, and 
glasses) over a wide range of state conditions 0, [1]. 
In particular, hard-sphere systems exhibit typical liquid- 
likephenomena, such as a first-order freezing transition 
[1, H @ and metastable glass states . 

Hard-hypersphere fluids (where the interaction poten- 
tial is infinite when two hyperspheres overlap and zero 
otherwise) are the natural extension of hard spheres to 
arbitrary dimensions d. Such systems have attracted an 
everlastin g at tention of many researchers 
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The main reason is twofold. First, studies 
of hard particles in high dimensions may reveal general 
behaviors of the equation of state (EOS), radius of con- 
vergence of the virial series, phase transitions, and fluid 
structure that can help to understand the corresponding 
properties in real fluids. Second, hypersphere systems 
provide well defined and very demanding test models for 
theoretical approximations to many-body problems. 

The thermodynamic and structural properties of d- 
sphere fluids in high dimensions have been examined by 
computer simulations %M,M.,M,M.,M,M,M.,M,d^- 
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At a theoretical level, a number of virial coefficients have 
been evaluated d, |T3, [H, [H, [H, HI], the asymptotic 
properties in the limit of inflnitely many dimensions have 
been investigated [l3jliB:j 22, 23, 26, 35|, several approxi- 
mate EOS have been proposed [l6|,|illi8j,(l9|,(23,l2l[2i, 
[27l . [28t , and scaled-particle and density- functional meth- 
ods have been applied to the fluid-solid phase transition 
[ll,!!^]- Regarding the structural properties, the Percus- 
Yevick (PY) closure to the Ornstein-Zernike (OZ) rela- 
tion has been proven to be exactly solvable at d = odd 
[1, [ll| , the solution having been worked out at d = 5 Q 
and d = 7 [3l|, IdSl (apart from the classical cases of d = 1 
and d = 3 [49l[50|). and overlap volume function rep- 
resentations have been proposed [Tl, [l3| • Additionally, 
special interest has focused on the packi ng pro blem and 
the formation of jammed structures [4ll. |43. |4^ l47l Isij . 



The primary aim of this paper is to introduce an an- 
alytical method for the study of hard-particle fluids in 
Euclidean spaces of arbitrary odd dimension. Specifl- 
cally, we derive an expression for a key function directly 
related to the static structure factor S{k) of the d-sphere 
fluid, from which all other structural and thermodynamic 
properties can be expressed. Our technique is based on 
the Rational Function Approximation (RFA) method, 
orig inally developed for three-dimensional hard spheres 
[52| and applied to a wide variety of problems [53|, in- 
cluding hard-sphere mixtures 15411 . sticky hard spheres 




square- well fluids [58!, l59|, and penetrable 
spheres 6£f\ . In the RFA approach developed in this pa- 
per we define a Laplace-space functional G{s) of the ra- 
dial distribution function g{r) that allows one to obtain 
S{k) in a simple way. By making use of the exact form 
of 5(7') to first order in density, a function ^'(s) directly 
related to G(s) is approximated by a Pade approximant, 
its coefficients being constrained by basic physical con- 
ditions arising from the small wavenumber behavior of 
S{k). The simplest implementation of the approach, i.e.. 
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the one with an equal number of coefficients and con- 
straints, turns out to coincide with the PY solution. The 
next extension contains two extra coefficients that are 
fitted to reproduce any desired EOS in a thermodynani- 
ically consistent way. 

The paper is organized as follows. Section |TT] summa- 
rizes the basic physical tools involved in this study and 
describes the application of reverse Bessel polynomials in 
the evaluation of Fourier transforms in odd dimensional 
space. An explicit expression for the overlap volume be- 
tween two hyperspheres, which plays a prominent role 
in the generalization of the RFA method, is derived in 
Sec, mil Section Hv] is devoted to the asymptotic expres- 
sions of S{k) and G{s) for long wavenumber, for short 
wavenumber, and for low densities. The generalization 
of the RFA approach to arbitrary odd dimension is pre- 
sented in Sec. IVland the corresponding evaluation of the 
direct correlation function c(r) is given in Sec. IVII Sec- 
tion IVIII shows explicit results for fluids in dimensions 
d < 11 and comparisons with available computer simula- 
tions for d — 5 and d = 7. The paper is closed with some 
concluding remarks in Sec. IVIIII 



II. FRAMEWORK 

A. Definitions 

The structure of a fluid is typically studied in terms 
of the radial distribution function, g{r), and a closely 
related function, the structure factor S{k), given by 



S{k) ^l + ph{k). 



(2.1) 



where 77 is the fraction of the total volume occupied by 
the c?-spheres (or packing fraction) and a the diameter of 
a particle. The general relation between p and 77 reads 



77 = Vdpcr 



Vd 



in/A) 



d/2 



r(i + d/2)' 



(2.5) 



where Vd is the volume of a d-dimensional sphere of unit 
diameter. In odd dimensions. 



Vd 



in/2) 



{d-l)/2 



dl\ 



(2.6) 



Furthermore, the structure factor is related to thermody- 
namics through the isothermal susceptibility as follows. 



(2J) 



The relation between the isothermal susceptibility and 
the compressibility factor is given by 



= -T-ivZ). 

dry 



(2.8) 



This equation can be used to impose thermodynamic con- 
sistency between the virial and compressibility routes to 
the EOS. 

For hard-particle systems in equilibrium, the structural 
properties considered here are athermal, i.e., the tem- 
perature does not play any relevant role. Moreover, the 
thermodynamic state can be characterized by a variable 
alone, i.e., the density, the pressure, or any related vari- 
able (e.g., 77, Z, x). 



Here k is the wavenumber, p is the density, and /i(fc) is 
the Fourier transform of the total correlation function. 



hir) = gir) - 1. 



(2.2) 



While gir) gives the relative probability of finding a par- 
ticle located a distance r from another particle located at 
the origin, S'(fc) is proportional to the scattered intensity 
of radiation from the fiuid and thus is obtainable from 
scattering experiments. An additional useful quantity to 
describe the fluid structure is the direct correlation func- 
tion c(r), which is deflned through the OZ relation. In 
Fourier space it reads 



c(fc) 



hik) 



1 



1 + phik) P 



1 



Sik) 



(2.3) 



The thermodynamics of hard d-sphere fluids can be 
fully accounted for by the compressibility factor Z = 
p/pksT (dimensionless combination of pressure p, den- 
sity p, and temperature T), which can be evaluated in 
turn from the contact value of the radial distribution 
fimction, g{a~^), 



Z=l + 2^"\gia+), 



(2.4) 



B. Low density expansions 

We consider briefly some exact low density results 
which will be included in the analytical theory presented 
in this paper. We start from the following general rela- 
tionship 



g(r) = [l + /(r)]y(r). 



(2.9) 



where 7/(r) is the cavity function and /(r) is the Mayer 
function, which in the case of hard d-spheres reads 



fir) 



-I, r < cr, 
0, r > a. 



(2.10) 



The low density behavior of g{r) can be derived from the 
virial expansion of the cavity function. 



yir) = 1 + yiir)p + y2ir)p^ + 



(2.11) 



where the functions yn(r) are represented by cluster di- 
agrams [ll, S]- In particular, the first-order contribution 
to the cavity function is 



7/1 (r)= /dr'/(r')/(|r-r'|). 



(2.12) 
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It is worth noting that, because of Eq. (j2.10|) . yi{r) 
represents the intersection volumen of two identical d- 
dimensional spheres of radius a whose centers are sepa- 
rated by a distance r [6ll |. 

With the preceding expressions, the lowest order terms 
in the 77-expansion of the radial distribution function, 

9ir)^9oir)+giir)7^ + 0{T^^), (2.13) 

are 

5o(r) = l + /(r) = e(r-a), (2.14) 



TABLE I: Reverse Bessel polynomials of degree less than 
seven. 

n Onjt) 

1 

1 1 + t 

2 3 + 3t + t^ 

3 15 + 15t + + 

4 105 + 105t + 45i^ + lOt^ + 

5 945 + 945t + 420*^ + 105t^ + 15t'^ + 

6 10395 + 10395t + 4725t^ + 1260f^ + 210t'' + 21*^^ + t'^ 



9i{r)^ [l + f{r)]y,{r), (2.15) 

where Q{x) is the step function [Q{x) = 1 if x > 1 and 
zero otherwise]. The factor 1 + /(r) prevents any pair of 
particle centers from getting closer than a distance a. 
Fourier transformation of Eq. (|2.15[) gives 



?i(fc) 



Vd 



yi{k) 



dk'yi(fc')/(|k-k'|) 



(2.16) 

where, by application of the convolution theorem in 



uk)= m 



(2.17) 



C. Fourier transform in odd dimensions 

The Fourier transform of an absolutely integrablc func- 
tion ^(r) in d dimensions is defined by 



(2.18) 



e(k) = J dre(r)e-^''-- 
and the associated inverse operation is given by 



1 



{2n) 



dk^(k)e*-. 



(2.19) 



where k is the wave vector. It is proven in Appendix 
Elthat, if the function ^(r) = £_{r) depends only on the 
magnitude r — |r| of the vector r and d = odd, then 
the d-dimensional direct and inverse Fourier transforms 
()2.18p and (|2.19p can be expressed as 

^(fc) = 2 ^ ld-2 Im{.F„K(r)](-»fc)} , (2.20) 



e(r) - 2 



(27r) 



-(d+l)/2 



-Im 



{^„K(fc)](-*r)}, (2.21) 



respectively. Here, Im(2:) denotes the imaginary part of 
z, n is an integer related to d by 



d — 3 

n = , a = 2n -|- 3, 



(2.22) 



and J^n[£.{x)]{s) is a functional of ^(x) defined by 

/"OO 

T.n[ax)]{s) EE / dxa;C(x)0„(sx)e-^^ (2.23) 



where the function 6n(t) is the so-called reverse Bessel 
polynomial of degree n [H, IH, Is^, IH, EH] , whose expres- 
sion is 



j=0 



(2n-j)! 
2"-J(rJ-j)!jr 



(2.24) 



Some of the properties of On{t) are summarized in Ap- 
pendix El Table H] provides the polynomials 0„(i) of de- 
gree less than seven. 

Two useful applications of Eq. (|2.23p correspond to 
£,{x) = 1 and to the step function £^{x) = 9(a; — 1). With 
the help of (|A5|) and (pTM]) one finds 



^n[l](.s)- 



0n+i{O) _ (2n+l)!! 



JP„[e(x-l)](s) = 



6'„+i(s)e" 



(2.25) 



(2.26) 



Henceforth we will indistinctly use d and n — (d — 
3)/2 in the remainder of the paper. Except for a few 
exceptions, we will generally follow the rule of employing 
n in subscripts and d in exponents. 



III. THE TWO-HYPERSPHERE OVERLAP 
VOLUME 

The intersection volume of two d-spheres whose cen- 
ters are separated a distance r is a key quantity in the 
study of hard systems [13, [67l| . Apart from its geometri- 
cal interest, it yields the cavity function to first order in 
density, yi{r), as mentioned below Eq. (|2.12p . Some ex- 
pressions of yi (r) in terms of special functions and recur- 
rence relations can be found in the literature [ij, [l^, [g^] , 
and an explicit expression has been recently derived by 
Torquato and Stillinger [42|. In this Section we use the 
representation (|2.2ip for d = 2n + 3 = odd to provide an 
alternative analytical expression of yi{r). 
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Taking into account that /(r) = 8(r — 1) — 1, where 
henceforth we take a as the length unit (i.e., a —X) and 
making use of Eqs. (H^D]), JMSl), and ^HI^ . one has [69|] 



(2^) 



(d-l)/2 



ik 



{ik)e 



-ikl 



(3.1) 

From the property (|A10[) it is easy to prove that, as ex- 
pected 

lim/(fc) = -^^^^^ = -2V (3.2) 



by 



The contribution yi{r) to the cavity function is given 



2/1 (0 



(2^) 



-(d+l)/2 



y.d-2 



dk kyi{k)9n(ikr)e 



'ikr 



where, according to (|2.17p and p.ip . 



J/2d 



[20n+l{ik)en+i{~ik) 



-2ik 



2ifcl 



(3.3) 



(3.4) 



Since j/i(0) = = finite, the integrand in (|3.3p 

is regular along the integration interval and so we can 
consider an integration path in the complex plane from 
k = — c» to fc = +00 that goes round the point fc = from 
below. The integral in Eq. (|3.3p decomposes into three 
contributions with integrands headed by e"*'^'', e~*'^'^''+^^ 
and ^-^Hr-'i) ^ respectively. If r > 2, we can close the 
path with a lower half-circle of infinite radius, so that the 
three contributions vanish. If < r < 2, however, the 
path in the third integral must be closed with an upper 
half-circle and the residue theorem yields a nonzero value. 
Therefore, for positive r one has 



= i?4„+4(r)e(2 



0, 



(3.5) 



where the residue 



R 



4n+4 



(r) 



Res 

fc=0 



fc-(2d-l)g|2_^^(„jfc)5)^(j/,^)gi/c(2-r-) 



is given by the term of order 2d — 2 ~ 4n + 4 in the 
fc-expansion of 9'^-^j^^{—ik)9n{ikr)e^^'^'^~'^\ i.e. 



R 



An+i 



(r) 



1 



(4n-(-4)! 



5t'"+'0^+i(-t)^„(rt)e(^-'-)* 



(3.7) 

where we have made the change t = ik. Equation (|3.7p 
implies that i?4„+4(r) is a polynomial of degree An -f 4, 
as indicated by the notation. Since 0^_|_]^(— t)0„(rt) is a 
polynomial of degree 3n + 2 in the Taylor expansion 
of the exponential factor contributes to Rin+i{'r) with 
factors (2 — ry with j > n + 2. Therefore, i?4„+4(r) 
factorizes into (2 — r)"+^ times a polynomial of degree 
3n -f 2. According to Eq. (|3.5p . the latter polynomial 



d-2 



^2ji+l 



starts with r' 
nite when r — > 0. From this analysis, we obtain 



because yi{r) must remain fi- 



Rin+iir) = 



(2-r)"+^P„+i(r), 



(3.8) 



where P„_|_i(r) is a polynomial of degree n+1. An explicit 
expression of Pn+i{r) is given in Appendix [Bl We finally 
obtain the first-order contribution to the cavity function 
for hard-sphere fluids in odd dimensions or, equivalently, 
the overlap volume of two hyperspheres of radius a = 1 
with centers separated a distance r (6l| , 



yi 



(r) =e(2-r)(2^)('^-i)/2(2_r)(<i+i)/2p^^^^(^), (3,9) 



Notice that yi{0) is equivalent to the volume of one hy- 
persphere the radius a- — 1, i.e., yi(0) = 2'^Vd- Moreover, 
setting r = 1 in Eq. (|B5[) yields 



yi(i) 
2/1 (0) 



1 



(2n + 3)!! 

2"+2 



(-4)- 



(2j + l)j!(n+l-j)!- 

(3.10) 

This result provides a simple expression for the third 
virial coefficient B3 of hard d-sphere systems in odd di- 
mensions. From the virial expansion 



1=2 



(3.11) 



and Eqs. l|2?4)) and (|2lT|l one has B2 = 2''"!?;^ and P3 
2'^-^Vdyi{l). Therefore, 



B3 ^o2/i(l) 
2/1 (0)' 



(3.12) 



IV. 



ASYMPTOTIC BEHAVIORS OF THE 
STRUCTURE FACTOR 



In this Section we examine the asymptotic long 
wavenumber, short wavenumber, and low density behav- 
iors of the structure factor and of a closely associated 
function, G(s), which will play a central role in the next 
Sections. We define G{s) as the functional Tn, intro- 
duced in (|2.23|) . of the radial distribution function. 



G(s) = Tn[9{r)]{s) - / drrg{r)0n{sr) 



(4.1) 



We have found that this is the optimal generalization 
to any odd dimension of the Laplace transform G(s) = 
£[rg{r)]{s) used in Rcf. [52"| in the study of hard spheres 
(d = 3, n = 0). With and (H^Sl), we note that 

G{s) = - 2-^-+/ drrh{r)en{sr)e-'''-. (4.2) 



Hence, using (|2.ip and (|2.20p . the structure factor is given 
by 

S{k)^l + p^^j^^i[G{ik)^G{~ik)]. (4.3) 
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A. Long and short wavenumber limits 

We shall now derive general conditions that the func- 
tions G{s) and S{k) must satisiy. First, one can eas- 
ily obtain an asymptotic expression of G(s) for large s 
by replacing rg{r) — 8(r — 1) [g(l^) + 0{r — 1)] and 
6»„(sr) = s"r" [l + 0{s-'^)] in (|tl]). The result is 



I.e. 



lim s^^-'^^^\'G{s) = 5(1+). 



(4.4) 



(4.5) 



It follows from Eqs. (|43)) and ((431) that, at long 
wavenumber, S{k) adopts the form 

^(^) « 1 + ;jSy7I^5(l+) sin [fc + J(5 - d)] . (4.6) 

Therefore, the structure factor for long k oscillates with 
an amplitude proportional to the contact value g{l~^), an 
envelope decaying as and a phase shift equal 

to 7r(5-d)/4. 

On the other hand, the asymptotic form of G{s) for 
small s can be derived from (|4.2p by using Eq. (|2.24p and 
the Taylor expansion of e~'^^. The result is 

(fj — 0]^ °° 

G{s) = + ^ a„,,i/,+lS^ (4.7) 



where 



dr h{r) 



(4.8) 



is the jVci moment of the total correlation function and 
the numerical coefficients a^j are given by 



E 



(4.9) 



One can verify that the first n coefficients a^j with j — 
2m + 1 = odd vanish, i.e., 



an,2m+l = 0, m = 0, 1, . . . ,71 - 1. 

Therefore, 

{d-2)U 



(4.10) 



G{s) 



2m 



values. Thus, inserting Eq. (|4.11|) into Eq. (|4.3p . we get 
the Taylor expansion of S{k) as 

oc 

5(fc) = 1 - 2'^d!!ry ^ (-l)™a„,2m+ii/2™+2fc2('"-"). 

m— n 

(4.12) 

Application of Eq. (|2.7p provides a direct relationship 
between the (d — l)th moment of the total correlation 
function and the isothermal susceptibility, namely 

X = 1 + 2''d7jHd-i, (4.13) 
where use has been made of the property 



(2n-H 1)!!' 



(4.14) 



B. Low density expansion 

The series expansion of the radial distribution function 
in terms of the packing fraction r?, Eq. (j2.13p . leads to a 
similar expansion for the hmction G(s), 



G(s) = Go(s)-KGi(s)ry + 0(772 



(4.15) 



with Gj{s) = J-'„[gj{r)]{s). The zeroth order term is 
derived from Eqs. (I?l4| and ((2?26)) . 



Go(s) = s-2g„+i(s)e- 



(4.16) 



Next, since, g{r) = 0(r — l)y(?') and yi(r) vanishes for 
r > 2, one can write 



1 r 

Gi(s) = — / drryi{r)0n{sr)e 

Vd J I 

It is proven in Appendix [Cl that 
A 



(4.17) 



= ^G^(S) + ^Q3n+4(S), (4.18) 



where 



A = (_i)('^-i)/22'^-id!! 



(4.19) 



and Q^n+iis) is a polynomial of degree 3n + 4 = (3rf — 
l)/2. Explicit expressions of the first few polynomials 
Q3n+A{s) are given in Table HIl 



V. RATIONAL FUNCTION APPROXIMATION 



+ «-2m+lff2m+2s2"+^ (4.11) 



The property (|4.10p is essential to guarantee that S{k) 
remains bounded at the limit of zero wavenumber and, 
consequently, the isothermal susceptibility takes finite 



The results presented in the preceding Sections are ex- 
act. In this Section, we propose the extension to hyper- 
spheres in arbitrary odd dimensions of the RFA, origi- 
nally introduced in the study of three-dimensional hard- 
sphere systems [H, [H, [T^I • The main steps in the RFA 
can be summarized as follows: (i) a functional G(s) of the 
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TABLE II: Polynomials Q3n+4,{s) associated with dimensions d = 3, 5, 7, and 9. 



n Q3n+4(s) 

12 + 12s - Gs" - 2s^ + (5/2)s* 

1 -2160 - 2160s - 360s2 + 360s^ + 305* - 42s'^ + 3s^ + (53/8)s^ 

2 1512000 + 1512000s + 4536008^ - 50400s^ - 47880s'' + 2520s'5 + 2940s^ - 300s^ - (285/2)s* 

+(789/16)s3 + (289/16)s"' 

3 -2667168000 - 2667168000s - 952560000s2 - 63504000s3 + 46720800s^ + 8618400s^ - 1738800s' 

-378000s'^ + 69930s** + 11130s3 - 3045s"^ - (585/128)s" + (38615/128)s^2 + (6413/128)s" 



radial distribution function g{r) is defined by a suitable 
Laplace transformation such that G{s) is simply related 
to the structure factor S{k); (ii) using as a guide the low 
density form of G(s), an auxiliary function ^{s) is intro- 
duced; (iii) the unknown function ^'(s) is approximated 
by a rational function (or Padc approximant), the degree 
difference between the numerator and denominator poly- 
nomials being dictated by the exact large-s behavior of 
G(s); (iv) finally, the coefhcients of the rational form for 
^{s) are determined by requiring consistency with the 
exact small-s behavior of G{s). 

According to Eq. (|4.3p . the first step described above 
is accomplished by the functional G{s) defined by Eq. 
(14. ip . Its large-s and small-s behaviors are given by Eqs. 
(14. 5p and (|4.4ip . respectively. In order to continue with 
the step (ii), let us rewrite Eq. (|4.I5p as 



„d-2 



„d-2 



G(s) 



[A-*i(s)e^]r,-f 0(77^), 



where, according to Eqs. (|4.16p and (|4.18p . 
^-0(5) 



(5.1) 



(5.2) 



(5.3) 



This suggests the introduction of the auxiliary function 
'^{s) through 



„d-2 



G{s) 



*(s)e'' - Ary' 



so that 



*(s) = «'o(s) + *l(s)77 + 0{7f). 



(5.4) 



(5.5) 



The large-s and small-s conditions (|4.5p and (|4.1ip imply 
that 

lim s-('^+i)/'*(s) = (5.6) 

s^oo g{l + ) 

'A77, j=0, 
0, l<J<d-l, 
di^{s)e^\^^^^{d\l{d-2)\\, j^d, (5.7) 
0, j ^d+l + 2m, 

TO = 0, . . . , n. 



respectively, where in Eq. (j5.7p we have taken into ac- 
count Eq. ((4T0)) . 

Thus far, Eqs. (|5.ip - (|5.7p are exact. Now we follow 
step (iii) and approximate the auxiliary function ^(s) by 
a rational form (or Pade approximant). 



*(s) 



BNjs) 
Am{s)'- 



M 

Am{s) = ^aj 



N 



s^ BN{s) = Y,bjs^. 
3=0 

(5.8) 

Note that the choice of a rational form for ^'(s) is com- 
patible with the exact requirements (|5.2p . (|5.3p . (|5.6p . 
and (|5.7p . The combination of Eqs. (j5.4p and (|5.8p con- 
stitute a simple approximation of G(s) that will be made 
consistent with the basic physical requirements outlined 
in Section IIVI To begin with, we note that the condi- 
tion (|5.6p fixes the relation between the degrees of the 
polynomials Am{s) and Bn(s), namely 



N = M ■ 



1 



(5.9) 



In fact, the ratio of the highest coefficients qm and 
directly gives the contact value of the pair distribution 
function: 



5(1^ 



bN 



(5.10) 



To close the RFA proposal (|5.8p we need to determine 
the coefficients {aj} and {bj}. Since one of them can be 
arbitrarily chosen, the number of independent unknowns 
is iV -I- Af -I- 1 . Following the step (iv), we resort to the 
exact small-s behavior (j5.7p , which imposes d -\- n -\- 2 = 
3n -|- 5 constraints. Therefore, in order to make the prob- 
lem solvable, one must have N + M > d + n + 1 = 3n + 4. 
In view of Eq. (|5.9p . this implies that 



N>d, M >n + l = 



1 



(5.11) 



The constraints stemming from Eq. (j5.7p or, equivalently, 
Eq. (|4.1ip are worked out in Appendix iDl 



A. Low density expansion 



It is worthwhile noting that the lower bounds (|5.1ip 
can also be derived by requiring consistency of the RFA 
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form ()5.8p with the exact zeroth order term, Eq. (|5.2|) . in 
a density expansion. In the approximation with the least 
number of unknowns, i.e., with 



N 



Eq. (15. 2p imphes that 



M 



bd 



r(=0 



fed 



1, 



(5.12) 



(5.13) 



i)=0 



where Sj^d is Kronecker's deha. As mentioned before, one 
of the coefficients in (|5.8p can be given a constant, arbi- 
trary non-zero value. Of course, this choice will not have 
any consequence on the physical results derived from 
(|5.4p . but an appropriate one may simplify the algebra 
involved in the further development below. In this sense, 
in view of (|5.13p . two adequate alternative choices are 
either 



Jn+i.,0 = (d-2)!! 



or 



(5.14) 



(5.15) 



Henceforth we will adopt (|5.14p . With the assumption 
(|5.12p and the choice (|5.14p . it is seen from Eq. (15. 



that the derivatives 



obey the following equations: 



^v^^'j in=o 



=0 and fo' = (9,,6j | 



»)=0 



(d-2)!! 



i-i 

m— max(0,j — M ) 

(5.16) 



ioT < j < N = d, and 



to the choice (|5.12p . We will refer to this case as the stan- 
dard RFA, i.e., the one in which all the Fade coefficients 
are determined from the basic constraints and so no free 
parameters remain. As will be shown in Sec. IVIi this 
standard RFA turns out to provide the exact solution of 
the FY integral equation for odd dimensions [8, 11] by a 
completely different method. 

For the case (|5.12p . Eq. (jPTp expresses the coefficients 
{bj} in terms of the coefficients {aj}. Insertion into Eqs. 
(|5.4p and (|5.8p allows one to rewrite G{s) in the form 



Gis) 



An+l{s) 



1 + >^vJ2'j=o aj(t)d-jis)' 
where we have called 



(5.18) 



E 



i-sy 



(5.19) 



Note that 0m (0) = 0. Equation (|5.18p assumes the nor- 
malization choice (I5.14P 71]. 

In order to evaluate the n + 1 coefficients {aj,j — 
1 , . . . , 71 + 1 } we can use the remaining constraints ob- 
tained from (jPip . An equivalent procedure results from 
the substitution of (|5.18p into (|4.7p . followed by a new 
power-series analysis. This latter method is simpler than 
the other one because the coefficients {bj} have already 
been eliminated. In this case, the analysis of powers 
yields 



1-2 



(5.20) 



where we have called 



M 



(5.17) 



for 1 < J < M = n-|-l, where qzn+i.j is the jth coefficient 
of the polynomial Q^n+iis)- Application of (|5.16p allows 
one to obtain recursively all the quantities 6^ , which can 
then be used in the evaluation of a' with the help of 

(Eni). 

In conclusion, the analytical form provided by Eqs. 
(|5.4p and (|5.8p exactly reproduces the zeroth- and first- 
order terms in density of G(s) when suitable values of 
the Fade coefficients and their derivatives are used. In 
the next subsection, we carry out the extension of this 
representation to arbitrary densities. 



B. Standard approximation 



Let us consider the RFA form (|5.8p when the number of 
unknowns (TV + Af + 1) equals the number of constraints 
(371 + 5). Taking into account Eq. ()5.9p . this corresponds 



(d-2)!!' 



7j 



(rf-2)!! 



(5.21) 



The coefficients Di are linear combinations of the {flj} 
given by Do = 1 and 



n+l 



A = yj - A77 ^ td-mjami I > 1, 



(5.22) 



where 



m=0 



_^ (-1)'+" 



(5.23) 



is the jth coefficient in the power series expansion of 
— e^0m(s). In Eq. (|5.20p and in the remainder of this Sec- 
tion we have adopted the conventions Oj = if j > 71 + 1 
and X^JLo • ■ ■ = if m < 0. 
Because of (|4.10l) . we have 



71 = 73 = • • • = 72n-l = 0. 



(5.24) 
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FIG. 1: Coefficients Uj, as a function of ry, obtained from the 
standard RFA metfiod for a fluid at d = 9. 



Therefore, Eq. (|5.20p with / — even <2n + 2 can be used 
to express the quantities 72™ with < m < n in terms 
of the coefficients Uj by means of the recursion relation 

m — 1 

72m = a2m+2 - -D2m+2 - ^ l2jD2(,n-j), < TO < n. 

3=0 

(5.25) 

The parameters 70, 72, 74, 76i ■ ■ - are linear, quadratic, 
cubic, quartic, ...in the coefficients {oj}. Next, Eq. 
with / = odd < 2n + 1 yields 



D 



2?ri+l 



2(m-j)-l — a2m+l, 



< TO < n. 



(5.26) 

When the {72m} obtained from Eq. (|5.25p are inserted 
into Eq. (I5.26P one gets a closed set of n + 1 equations 
for oi, 02, . . . , fln+i- Therefore, the implementation of 
the standard RFA method reduces to solving a set of 
n + 1 = (d — l)/2 algebraic equations, which become 
nonlinear for n > 1 or, equivalently, d > 5. In general, 
the number of mathematical solutions (including complex 
ones) is 2". In the case of multiple solutions, we choose 
the solution which yields the correct low density limit 
given by (|5.13p . In this sense, it is important to note 
that the asymptotic behaviors of aj and bj as 77 — > 
[Eqs. (|5?T3l) . jSHH), and (fSTZ]) ] are naturally included 
among the solutions of Eqs. (|5.26p . This is due to two 
reasons: (a) Eq. (|4.1ip is verified by the exact density 
expansion of G(s) up to any order, in particular to zeroth 
and first orders [Eq. (|4.15p ]: (b) the number of constraints 
derived from (|4.1ip coincides with the least number of 
Pade coefficients required in the exact expansion of G(s) 
up to first order in density. 



The set of equations (|5.26p can be analytically solved 
for d = 1, 3, 5, and 7, the solutions for d = 1, 3, and 5 
being explicitly given in Appendix lEl For d = odd > 9, 
however, the solutions must be obtained numerically. By 
means of explicit evaluations in all the cases analyzed, 
we have found that the coefficients {aj} are finite for 
< ?7 < 1 and adopt the following values at rj — 1: 



2J ■ 



(5.27) 



Figure [T] illustrates the physical roots of Eq. (|5.26p for 
a fluid of dimension d = 9 (n = 3) as functions of the 
packing fraction. It may be observed that, at the low 
density limit (77 — > 0) , each aj tends to the corresponding 
coefficient uj^j of the reverse Bessel polynomial 64(1) (see 
Table m, as dictated by Eq. (|5.13p . These roots [except 
ao which is fixed by Eq. ()5.14p ] monotonically decrease 
as the density increases, reaching their minimum values 
given by Eq. (|5.27p at 77 = 1. 

Once the coefficients {aj} are determined as functions 
of 77, the structure factor of the system is given by Eqs. 
(14. 3p and (|5.18p . The ansatz (|5.8p thus provides an ex- 
plicit expression for the structure factor to any finite den- 
sity. Because the approximation to G(s) is exact to first 
order in density, the structure factor thus obtained is ex- 
act to second order in density. In conclusion, the 'in + 5 
constraints derived from the small-s behavior of G(s), to- 
gether with the zero-density conditions (|5.13p . allow us 
to completely characterize S{k) for a hypersphere fluid 
in odd dimension in the standard RFA approach. Since 
S{k) and h{k) are directly related through Eq. (|2.ip . ap- 
plication of Eq. (|2.2ip yields the radial distribution func- 
tion 5(7-). 

The compressibility factor resulting from the virial 
EOS (|2.4p is determined by the standard RFA method 
through the relations (|5.10p and (jD7p . The result is 



Z.,{f]) = 1 + 2''-i7/a„+ 



1 + A77E^"'^ 



id-j) 



(5.28) 

The isothermal susceptibility x given by (|4.13p can also 
be easily evaluated. On the one hand, from Eq. (|4.14p 
and the definition of the 7j factors [Eq. (|5.2ip ] one has 

i7,_i = (-l)«+i[(d-2)!!]27,_2. (5.29) 

On the other hand, Eq. (|5.20p at / = d yields 



7d-2 



Dd +y^72j£'rf-2-2j 
3=0 



(5.30) 



Therefore, from Eq. (|4.13p we obtain 

n 

Xiv) = 1 - 2A(d - 2)!!7y(i?d + ^ 72,i?d_2-2,) , (5.31) 



j=o 
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with factors 70,72, ■• ■ ,7d-3 given by (15.2511 . Once x is 
known, the thermodynamic relation (12. 8p can be inte- 
grated to obtain Z in the so-cahed compressibihty route. 



(5.32) 



Additionally, one finds from (I5.27P the following limit 
values 



Z,(l) = c^, x(l)=0. 



(5.33) 



Comparison between the results obtained from Eqs. 
(|5.28p and (|5.32p . which gives a measure of the degree of 
thermodynamic inconsistency in the standard RFA solu- 
tion, will be presented in Sec. I VIII 



C. Extended approximation 

It is possible to construct RFA solutions more elabo- 
rate than the standard RFA one by considering in Eq. 
(15. 8p a number of unknowns iV + A/ + 1 larger than the 
number 3n -|- 4 of basic constraints (|5.7p . We will refer 
to this case as the extended RFA. The simplest extension 
corresponds to 



A^ = d+1, M = n 



(5.34) 



since Eq. (j5.9p must be preserved. This involves two new 
parameters (a„+2 and hd+i) which can be freely chosen 
without compromising the basic physical requirements. 
A natural choice is to adjust a„+2 and hd+i by requir- 
ing prescribed values of the contact value of the radial 
distribution function, gc = g{l^), and of the isothermal 
susceptibility x- In practice, one can use p.Sp to eval- 
uate X from gc (or vice versa), so that only one EOS is 
needed and the thermodynamic consistency between the 
virial and the compressibility routes is thus guaranteed 
by construction. 

With the choice (|5.34p . the relationship (|5.10p can be 
used to eliminate one of the two new Fade coefficients (for 
instance, hd+i) in terms of the other one (a„+2). Making 
use again of Eq. (|D7p . the extended RFA for G(s) reads 



G{s) 



1 + gc ^an+2S + A77 Yn=o a,j(l)d~j{s) 



(5.35) 



3=0 



The set of Fade coefficients , < j < n + 2}, is ob- 
tained in a manner similar to that given above. Thus, oq 
remains fixed by (|5.14p . while ai, . . . , a„+i are related 
to an+2 by the 71 -I- 1 equations (|5.26p . So far, by setting 
o-n+2 = we recover the standard RFA. However, now 
we fix the prescribed x and then Eq. (j5.3ip provides the 
needed equation to close the set (|5.26p . 

The nonlinearity of the problem in the extended RFA 
is higher than in the standard case, the number of math- 
ematical solutions being 2"+^. In particular, one has to 
deal with a quadratic equation for d — 3 [52, [El, [lO] , a 



quartic equation for d = 5, and so on. It can be verified 
that in the zero-density limit the physical solution has 
the form 



"'j\r,=0 — ^n+lj + an+2\rf^Qt^n+l,j-l, 



(5.36) 



^jlr,=0 = ^3,d + 0^n+2\r,^oSj.d+l, (5.37) 

where the numerical value of an+2|^^o depends on the 
value of the fourth virial coefficient predicted by the pre- 
scribed EOS. Equations (|OS)l and ((07)l imply that 

lim 74„+2(s) = 9n+i{s) (^1 + a„+2|^^o^) ' (5.38) 



lim Bd+iis) = 
77— >o 



1 + an+2 



1^=0 . 



(5.39) 



so that Eq. (|5.2p is recovered, irrespective of the chosen 
EOS (provided, of course, it is consistent with the exact 
second and third virial coefficients). 

It is important to note that the physical root of the 
set of equations (|5.26p and (|5.3ip must correspond to 
an+2 > 0. Otherwise, since gq — (d — 2)!! > 0, there 
would exist at least one positive real root sq of the poly- 
nomial A„_|_2(s). According to Eq. (|5.35p . this would 
imply G{so) — 0, what is incompatible with a positive 
definite g{r). A careful analysis of Eqs. (|5.26p and (|5.3ip 
shows that the condition an+2 > requires that the cho- 
sen values of Z and x satisfy the inequalities Z > Z^^ 
and X > Xc^7 where Z^"^ and are the compressibil- 
ity factor and the isothermal susceptibility, respectively, 
obtained from the standard RFA or, equivalently, from 
the FY solution. Therefore, if the prescribed x is ob- 
tained from the prescribed Z by application of Eq. (|2.8p , 
the extended RFA provides physical correlation functions 
only if 



7PY 



<Z <Z. 



FY 



(5.40) 



When this condition is verified, the physical root of the 
set of equations (|5.26p and (|5.3ip corresponds to the 
smallest positive real value of an+2- 



VI. DIRECT CORRELATION FUNCTION 



According to Eq. (jA12p , the direct correlation function 
c(r) can be evaluated from its Fourier transform as 

c(r) = ■5-^4=2 M dkkc{k)0n{ikr)e-'''''. (6.1) 



In the RFA approach, using Eqs. (gSl), (|5^ . and ([O)) . 
we obtain 



S{k) 



B+B- - A+A- 



(5+e«fc - A+)(B-e-»'= - A-)' 



(6.2) 
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where = XrjAMiHk) and = BN{±ik). Next, use 
of Eq. ([231) yields 



Since, according to Eqs. (|4.7p and (|5.4p . one has 
BN{s)e'' - A?7Aa/(s) = s'* [1 + 0(s)], it follows that the 
denominator in Eq. (|6.2p is of order fc^''. Therefore, the 
numerator must also be of order k^"^ to have a finite value 
of S{0). More specifically, for both the standard {N = d, 
M = n+1) and the extended {N ^ d + 1, M = n + 2) 
RFA approaches, one has 



B+B- 



A+A- 



{bd-ib. 



'd+l 



2d 



(6.4) 



A. Standard RFA 

In the case of the standard RFA, one has hd+i = 0, so 
that B+B^ - A+A^ = -bl{ik)'^'^. Consequently, when 
inserting Eq. (|6.3p into Eq. (|6.ip . each one of the three 
integrands has a pole of order 2d — 1 at fc = 0. A residue 
analysis similar to that employed in Sec. IIIII leads to 

cW = i?4„+4(r)e(l - r), (6.5) 



where R4n+i{r) is a polynomial of degree 4n + 4 : 
given by 



i?4n+4 (r) = Res 



An+i{-s)Bdis) 



i(sr) 



2d-2 



(6.6) 



Comparison with Eq. (|3.5p shows that, in the standard 
RFA, c(r) has in the region < r < 1 a polynomial 
form similar to that of yi{r) in the region < r < 2. 
As happened in the latter case, finiteness of c(0) implies 
that c(r) is a polynomial of degree 4n + 4 — (d — 2) = d. 

Equation (16. 5p shows that c(r) = for r > 1 in the 
standard RFA. Of course, the standard RFA complies 
with the physical requirement g{r) =0 for r < 1. These 
two conditions define precisely the FY closure to solve 
the OZ equation for hard hyperspheres. Therefore, we 
find that the standard RFA, i.e., the approximation given 
by Eqs. (|5.4p and (|5.8p with the least number of coeffi- 
cients satisfying the requirements (|4.5p and (|4.1ip . coin- 
cides with the FY solution for hard hyperspheres of odd 
dimensions |Tl| . This is a remarkable result since both 
approaches are in principle rather independent. In fact, 
following the philosophy behind the RFA approach, it is 
straightforward to proceed to the first natural extension 
of the FY solution or standard RFA. 



B. Extended RFA 

In the extended RFA approach, bd+i ~ 7^ 
0. As a consequence, Eq. ()6.4p shows that, besides the 



pole of order 2d — I at k = 0, each one of the three 
contributions to kc{k) has also two simple poles k = ztiK, 
K being a real quantity given by 



Therefore, applying the residue theorem we get 
c{r) = i-ir+^y\K^ + K_ Ko)'-^^^^ 



(6.7) 
(6.8) 



for r > 1 and 



c(r) 



(_l)(d+l)/2 
„d-2 



for r < 1. In the above expressions, 

^ An+2iTn)Bd+i{±K)e^'' 



_ , An+2{K)An+2{-li) 



°d+l'^ 



2d 



(6.9) 

(6.10) 
(6.11) 



^4n+4('^) 



Res 

s=0 



An+2{-s)Bd+l{s) 



o2d- 



l(s2 



i{sr)e 



(6.12) 

Equation (|6.8p shows that, in contrast to the FY solu- 
tion (or, cquivalently, the standard RFA), c(r) does not 
vanish outside the core in the extended RFA. In fact, 
c(r) oc 6'„(Kr)e~'"'/r'^~2 for r > 1. This functional form 
can be considered as an extension to odd d of the well- 
known Yukawa form e~'^"^ jr for d — i [73|. 

The form (|6.8p obtained from the extended RFA for 
hard hyperspheres admits an alternative interpretation. 
Imagine an interaction potential with a hard core at r = 1 
plus an attractive Yukawa tail —eOn{zr)e~'^'^/r'^~'^ for 
r > 1. In the mean spherical approximation (MSA) 
P] the closure to the OZ relation would be c(r) = 
{e/kBT)9n{zr)e~^^ /r'^^'^ for r > 1, which has the same 
form as Eq. (|6.8p . As a consequence, the MSA for the 
hard-core Yukawa interaction with d = odd is exactly 
solvable, the corresponding functions G{s), S{k), and 
c(r) for r < 1 being given by Eqs. (|5.35p . (|6.2p (with 
N = d + 1 and M = n + 2), and (p^ . respectively. 
The main difference between the extended RFA for hard 
hyperspheres and the MSA for the hard-core Yukawa po- 
tential is that, while in the former case the two input 
parameters are the contact value gc — 3(1^) and the 
isothermal susceptibility x, in the latter case the control 
parameters are the reduced temperature ksT /e and the 
inverse interaction range z. 

In the three-dimensional case, the extended RFA [s^l 
reproduces the the so-called generalized mean spherical 
approximation (GMS A) for hard spheres [t^ , which con- 
sists of closing the OZ equation with the assumption that 
c(r) has a Yukawa tail outside the core (r > 1). From 
that point of view, the extended RFA applied to higher 
odd dimensions can be seen as the natural extension of 
the GMS A to those dimensions. 
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Packing fraction, -q Packing fraction, 77 



FIG. 2: Top panel: Compressibility factors and 
(curves on the right and left edges of each shaded area, respec- 
tively) arising from the virial and the compressibility routes 
according to the standard RFA method (PY theory) at di- 
mensions d — 1, 3, 5, 7, 9, and 11, as functions of the packing 
fraction ry. Bottom panel: Difference between Z^^ and Z^^ . 



FIG. 3: (Color online) Compressibility factors Zc and Zv 
arising from the standard RFA (edges of shaded areas) for 
d = 5, 7, 9, 11. The values obtained from the SMS EOS ^ 
for d = 5, 7, 9, 11 (triangles), and the LM EOS l9] for d = 5, 7 
(squares) are also shown. The vertical dashed lines mark the 
freezing packing fraction reported in Table Hill 



TABLE III: Packing fractions at freezing transition {rff) and 
closest packing (rjc) predicted for fluids at d = 1, 3, 5, 7, 9, 11. 



d 



1 
3 
5 
7 
9 
11 



0.494" 
0.19'' 
0.072" 
0.027^* 
0.009^* 



V2n/6 ^ 0.74048" 
v^ttVSO ^ 0.46526" 
7rVl05 ^ 0.29530" 
^7r''/945 « 0.14578" 
327rV(93555\/3) ^ 0.06043" 




"Monte Carlo simulations 0] 
''Molecular dynamics simulations 
"Molecular dynamics simulations 
"^Estimated by the method of Refs 
"^Densest lattice packing listed in Ref. 



VII. RESULTS 

A. Standard RFA (PY theory) 

We list in Appendix [E] the explicit solutions of Eq. 
(|5.26[) corresponding to the standard RFA approach (or, 
equivalently, the PY theory) for fluids in dimensions d ~ 
1,3,5. The solution for d = 7 is also analytical and 
agrees with the results reported in Refs. [s!], |4^. The 
solutions to Eq. (|5.26p for d = 9, 11 have been obtained 
numerically. 

The compressibility factors and derived by 
the standard RFA approach from the virial and com- 
pressibility routes [Eqs. (|5.28|1 and (j5.32|) ]. respectively. 



are shown in Fig. [2]for d = 1, 3, 5, 7, 9, 11, as functions of 
the packing fraction. Both routes yield identical and ex- 
act results only in the case d=\. Discrepancies between 
and grow noticeably with increasing dimension 
for d > 3 (bottom panel in Fig. (2]). The compressibility 
factors predicted by the standard RFA have a singularity 
at 77 = 1 for all d [cf. Eq. (|5.33p ]. However, since d-spheres 
are not space filling (except for c? = 1) the true pressure 
must present a singularity at a certain density lower than 
or equal to the closest packing fraction 77c. The maximal 
packing fractions presently known in these dimensions 
are listed in the third column of Table IIIII In practice, 
the inconsistency between Z^^ and Z^^ is not as severe 
as it is appears in Fig. [2] because the fluid phase is actu- 
ally limited to very low values of 77 for high dimensions. 
Freezing transitions are predicted for hard d-sphere fluids 
with d> \ and found to occur at increasing lo wer p acking 
fractions with increasing dimensionality [6[[30l.l4l|. Some 
values of the freezing packing fraction r\f are listed in the 
second column of Table IIIII Figure [3] depicts the com- 
pressibility factors Z^^ and Z^^ up to densities in the 
neighborhood of the liquid-solid phase transition for flu- 
ids with d = 5, 7, 9, 11. The figure also includes the values 
of Z predicted by the extension of the Carnahan-Starling 
EOS proposed by Song, Mason, and Stratt (SMS) [l^, 
as well as the predictions of the semi-empirical EOS pro- 
posed by Luban and Michels (LM) [l^. There is a rea- 
sonable agreement among these two latter EOS and the 
PY results via the compressibility route within the fluid 
phase. In fact, comparison of Z^^ and Z^^ with com- 
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10 15 5 10 15 

wave number, k 

FIG. 4: (Color online) The structure factor as obtained from 
the standard RFA (or PY solution) for hyperspheres in d = 5 
(at ri = 0.06,0.13,0.2), d = 7 (at = 0.03,0.05,0.07), 
d = 9 (at 77 = 0.01,0.025,0.04), and d = 11 (at 77 = 
0.005,0.01,0.015). The first peak increases with increasing 
density. 




1 

0.5 


-0.5 



d = 5 



I ' ' ' ' I 
d = 7 
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FIG. 5: (Color online) Top panels: Direct correlation function 
for d = 5 (at ry = 0.2) and d = 7 (at rj = 0.07) as computed 
from the standard RFA (solid lines) and the extended RFA 
using the SMS EOS [l3| (dotted lines) and the LM EOS ^ 
(dashed lines). Bottom panels: Details of the Yukawa tails 
outside the core as obtained from the extended RFA solutions. 



puter simulations in the cases d = 5 [1^ and d = 7 (3]| 
shows that, in the stable fluid region, the true values of 
Z are bracketed by Z^^ and Zf^ in the form indicated 
by Eq. dHOl), being closer to Zf^ than to Z^^ . This 
in passing gives further support to the extended RFA, 
which requires as input an EOS satisfying Eq. (|5.40[) in 
order to provide physically meaningful results. 

Figure [J] shows the structure factor obtained from 
the standard RFA method for fluids in dimensions d = 
5, 7, 9, 11 for densities near the values of 77/ given in Table 
mil As expected, the oscillations of S{k) become more 
pronounced as the density increases. We can also observe 
that the location of the first peak tends to move to higher 
wavenumbers as the dimensionality increases (4l| . 



B. Extended RFA 

As seen in Sec. IVIl the main signature of the extended 
RFA solution with respect to the standard one (or PY 
solution) is that the former yields a direct correlation 
function with a (generalized) Yukawa form outside the 
core [cf. Eq. (|6.8p ]. In Fig. [5] we compare c(r) evaluated 
from the standard and extended RFA (the latter being 
complemented by the SMS and LM EOS) for d = 5 and 
d — 7 at densities close to the respective values of 77/. 
At those densities, as shown in Fig. [31 the SMS EOS 
is closer to the PY compressibility route than the LM 
EOS. As a consequence, the value of the extended co- 
efficient a„_|_2 is smaller in the former case than in the 
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radial distance, r 

FIG. 6: (Color online) Radial distribution function of hard hy- 
perspheres in d = 5 obtained from the extended RFA method 
(dashed line), the PY solution (solid line), and Monte Carlo 
simulations ^3% (symbols) at p = 0.8 (7; ~ 0.1316). The 
top panel shows the global behavior, while the bottom panels 
show the details near the first peak (left) and between the 
second and third peaks (right). 
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FIG. 7: (Color online) Radial distribution function of hard 
hyperspheres in d = 7 obtained from the extended RFA 
method (dashed line), the PY solution (solid line), and molec- 
ular dynamics simulations [iol . [4^ (symbols), at p = 1.4 
(ri ~ 0.05168). The top panel shows the global behavior, 
while the bottom panels show the details near the first peak 
(left) and near the first minimum (right). 



latter. This explains the fact that in Fig. [5] the curves 
corresponding to the extended RFA complemented by the 
SMS EOS are much closer to the PY ones (and with a 
weaker Yukawa tail) than those obtained by using the LM 
EOS. Given the semi-empirical character of the LM EOS 
and its excellent agreement with computer simulations 
[igl . [3l| , it is reasonable to expect that the true c(r) is 
better represented by the extended RFA complemented 
with the LM EOS than with the SMS EOS. Of course, 
there are other EOS proposed in the literature for hard 
hyperspheres [Hj [ICi, .25], but the ones considered in 
Fig. [5] are enough for illustrative purposes. 

Comparisons between the pair distribution function 
obtained from both RFA approaches and available com- 
puter simulations [13, SO, IS] are shown in Figs. [H] and 
[7] for d = 5 and d = 7, respectively. The extended RFA 
results have been obtained with the SMS EOS, but we 
have checked that no significant differences are found if 
the LM EOS is used instead. This indicates that g{r) is 
much less sensitive than c(r) to the choice of the input 
EOS. Figures [6] and [7] show that the standard RFA (or 
PY solution) provides an accurate estimate of g{r) for all 
radial distances at the considered densities (77 w O.777/), 
although some small discrepancies appear near the first 
maxima and minima, especially at contact. These devi- 
ations are substantially corrected by the extended RFA. 



VIII. CONCLUDING REMARKS 

In this work we have generalized the RFA method, 
originally developed for three-dimensional hard-sphere 
fluids [s^, m, [201 1 *o hypersphere systems in arbitrary 
odd dimensions d, providing explicit results for d < 11. 
This generalization is not trivial at all. In the application 
of the RFA approach, one must define a functional G{s) 
of g{r) in Laplace space, introduce an auxiliary function 
^1^(5) directly related to G{s), and approximate ^'(s) by 
a rational function (or Pade approximant), determining 
the coefficients by the application of basic consistency 
conditions. 

In the one- and three-dimensional cases the key func- 
tion G(s) is the Laplace transform of g{r) and rg{r), 
respectively [H, HI] . Thus one might be tempted to de- 
fine G{s) for d = odd simply as the Laplace transform of 
7.('i-i)/2^(^)^ Here, however, we have adopted the crite- 
rion that G(s) must be defined as to be closely related 
to the static structure factor S{k) of the fluid. This has 
led us to the definition (|4.ip and to the relationship (|4.3p . 
where the reverse Bessel polynomials (|2.24p play a central 
role. As a byproduct, we have derived a general polyno- 
mial expression for the overlap volume of two identical 
hyperspheres, this quantity providing the cavity function 
to first order in density, yi{r). Once G{s) has been iden- 
tified, one needs to define the auxiliary function ^(s) to 
be approximated by a rational function. Using the exact 
knowledge of G{s) to first order in density, Eqs. (|5.1[) - 
(15. 3p . it turns out that the natural definition of 4'(s) is 
provided by Eq. (j5.4p . Finally, ^'(s) has been approxi- 
mated by a rational function, Eq. (|5.8p . the degree dif- 
ference between the numerator and denominator being 
fixed by the exact large-s behavior of G{s). The coeffi- 
cients in the Pade approximant for ^'(s) are constrained 
to fit the exact small-s behavior of G(s) or, equivalently, 
the small-fc behavior of S'(fc). This representation repro- 
duces the exact expansion in density of S{k) up to second 
order. 

We have called standard RFA to the case in which the 
number of parameters in the Pade approximant equals 
the number of constraints. It turns out that the associ- 
ated direct correlation function vanishes outside the core. 
Therefore, quite interestingly, the standard RFA coin- 
cides with the solution of the PY closure to the OZ rela- 
tion for hard-particle fluids in all odd dimensions d, [Tl| . 
This equivalence between two completely independent 
paths allows one to view the PY solution for hard hyper- 
spheres as the simplest one of a broad class of approxima- 
tions. In fact, a more flexible approximation is obtained 
by adding a pair of new terms (one in the numerator 
and the other one in the denominator) in the Pade ap- 
proximant for ^'(s), resulting in what we have called the 
extended RFA. Apart from satisfying the small-fc behav- 
ior of 5(fc), the parameters are determined by requir- 
ing thermodynamic consistency with a prescribed EOS, 
which must satisfy the inequalities (|5.40[) to ensure that 
the extended RFA solution is positive definite. 
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Comparison with available computer simulations for 
d = 5 [sT] and d = 7 [13, IZ^I shows that the radial 
distribution function predicted by the standard RFA (or 
PY solution) is rather accurate. On the other hand, there 
exist certain small discrepancies (especially near contact) 
that are satisfactorily corrected by the extended RFA. 

The work presented in this paper is aimed at contribut- 
ing to our understanding of the structural properties of 
hard-hypersphere fluids and to the mathematical intri- 
cacies of their statistical-mechanical description. More- 
over, this work paves the path to the study of other re- 
lated systems in d dimensions, such as sticky hard hyper- 
spheres, square-well particles, or multicomponent hard- 
hypersphere fluids. Work is now in progress along these 
lines and the results will be published elsewhere. 



as follows: 



jn(x) 



-ix)e^^ — 9n{ix)e 



— Im[0„(-zx)e"], 



(A4) 



where the reverse Bessel polynomial 0„(t) is defined by 
Eq. (12:241) [765. 

A useful integral identity for On(t), which we have not 
found in the literature, is 



On{t)e-'^ J dzz0„_i(z)e-^ (A5) 



By taking the derivative of both sides of Eq. (jASP one 
obtains 
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APPENDIX A: FOURIER TRANSFORM USING 
REVERSE BESSEL POLYNOMIALS 



If the function ^(r) — ^{r) is isotropic, Eqs. (|2.18p and 
(Pl^ become [75] 



(Al) 



(27r)'^/2 7o 



respectively. Here k — |k| is the magnitude of the wave 
vector and J,y{x) is the Bessel function of the first kind of 
order v. For half-integer order v = n + 1/2, it is usual to 
introduce the spherical Bessel function of the first kind 
jn{x) given by 



Next, using Eq. (jA4p . it is easy to prove that the recur- 
rence relation of the spherical Bessel functions [ttI ]. 

4-jnix) = -^^^-^jn{x) + jn-lix), (A7) 

da; X 



is recovered. In passing, from Eq. (|A6p we may note the 
recurrence formula 



dt 



(A8) 



Another recurrence relation is [gJI 

0„(O - (2n- l)0„_i(i) +i2(?„_2(t). (A9) 



Besides, with the help of Eq. (IA4I) and expression (10.1.2) 
of Ref . 77] , one can find the asymptotic relation 

„■ ^2n+l 

- [9^{zx)e-^^ - 9,,{-^x)e^'■^] = [l + 0{x^)] , 

(AlO) 

which is used in the evaluation of yi{r) in Sec. IIIII 

We return now to the problem of expressing Fourier 
transforms in odd dimensions. The functions ^(r) and 
£,{k) can be extended to negative r and fc as ^(— ?■) = ^(r) 
and ^(— fc) ~ £,{k), respectively. Then, with the help of 
(|X3)) and (M^ we can rewrite and ((X2|) as 

g(fc) = ^ » / drrC(r)0„(^fcr)e-^'=^ (All) 



g(r) = ^ 1 I AkkmOu{ikr)e-'''\ (A12) 



in(x) — \/ ^^,1+1/2(2;). 



(A3) 



In order to rewrite the Fourier transform in a more ap- 
propriate form, we express the spherical Bessel functions 



respectively, where n is defined by Eq. (|2.22p . Finally, in- 
troducing the functional (|2.23p . one arrives at Eqs. (|2.20p 
and dMI]). 

It is worthwhile noting that if 9n{x) is replaced by 
its polynomial expression (I2.24p . then !Fn[£,{x)] can be 
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expressed in terms of the Laplace transforms of x^^^^{x) 
with < J < n, 

n 

J'nia^m = ^C.„,,,S^£[.T^"+1C(X)](S), (A13) 

3=0 



m^m^ / dxax)e-'^. (A14) 

JO 



APPENDIX B: EVALUATION OF P„+i(r) 

In this Appendix we obtain an exphcit expression for 
the polynomial P„+i(r) related by Eq. (|3.9p to the first- 
order cavity function yi{r). The derivation of P„+i(r) 
from relations (|3.7|) and p.8|) is straightforward but 
rather tedious. We found it more practical to use a known 
expression for the scaled overlap volume in terms of the 
normalized incomplete beta function 16], 



yi(0) S((d + l)/2,l/2) 



where 



B^{a,b) / dtt''~'^{l - t) 



\b-l 



is the incomplete beta function and B{a, b) ~ 
Bx=i{a,b) = r(a)r(6)/r(a + b) is the beta function. In 
the present case, for d odd, 

B{{d + l)/2, 1/2) = 2^^+^^/^((^-l)/2)! ^ ^^3^ 



From (|Bl[) " (jB3p . it is straightforward to obtain 



dr 



^2X (d-l)/2 

2(<i+i)/2((rf _ l)/2)! " T j e(2-r). 

(B4) 

Expanding the binomial and integrating over r, one has 

"+i (-rV4p 



Qdir) = 1- 



(2n + 3)!!r 



2n+2 (2j- + i)^-!(„ + i_j)r 



(B5) 



where it has been implicitly assumed that r < 2. 

On the other hand, inserting Eq. p.9p into the defini- 
tion of ^d{r) we get 

n,ir) = (^!i±^(2 - r)"+2p„+i(r), (B6) 

where we have taken into account that, since yi{0) = 
2%, thenP„+i(0) = 2-("+iV(2n + 3)!!. From Eqs. jHI 
and (jB6[) one finally gets 



n+l 



j=0 



(B7) 



with the coefficients 



2-(m+j) 



j!(2TO + 1)!! 



[(i-i)/2] 

E 

(=0 



(-l)'(m-l+j-20! 
(2; + l)(j-l-20!Z!(TO-0! 



(B8) 

where [(j — l)/2] represents the integer part of {j — l)/2. 

Equation (|B5P coincides with the expression (3-24) of 
Ref. [l^l evaluated in the case d — odd. The novel con- 
tribution of the procedure outlined in Sec. Illll is to show 
the factorization of the overlap volume [Eqs. p.9p and 
(jB6P ] into the product of (2 — r)"+^ and P„+i(r), which 
is not evident from Eqs. ((BT|) . (lB4l) . or (|B5|) . 



APPENDIX C: EVALUATION OF Gi(s) 

In this Appendix the function G'i(s) defined by Eq. 
()4.17p is evaluated. We start from the identity 



drr"(sr)^e-- = -^(m+j)!^ 



(sr)' 



(CI) 



1=0 



(B2) Consequently, 



drr™^?„(sr)e-^'- = ^ (m + j)!c^„ 



X E 



Finally, making use of Eq. (|B6p . we get 

22n-H3 



drryi(r)6l„(sr)e '''' = -^5^Q3n+4(r, s), 



2/i(0) 

where we have called 

g3„+4(r,s) EE (2s)2"+-^5](j + l)!c.„,,X^ 

n+l 



(C3) 



(2n-H3)!!:^' (-I)"(4s2)»+i- 



(2to -|- l)TO!(n -|- 1 — m)! 



n 2m+2+j , w 

x5^(2m + 2 + j)!c.„,, ^^C4) 

Equations (j4.17p and (|C3P allow us to write 

Gi(s) = s-(2"+5) [e--Q3„^4(i,s) _ e-2-g3«+4(2,s)] . 

(C5) 

The function (53,1+4(2, s) can be further simplified. To 
that end, note that, according to Eq. (IC3p . 

^Mr)\{s) = [g3«+4(0, s) ~ e-2-g3n+4(2, s)] . 

(C6) 
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Using Eq. (|2.20|) and comparing with Eq. (|3.4p one ob- 
tains 

Q3n+4(2, s) = (-l)"4"+i(2n + mOl+As), (C7) 



?3n+4(0,s) + Q3n+4(0,-s) = (- 1)"22"+3 (2,1 + 3) !! 

x6l„+i(s)e'„+i(-s). 

(C8) 



Equation (|4.18p is obtained by inserting Eq. (jC7|) 
into Eq. jCU, using Eqs. (|4?T6ll and (|4?T9l) . and calling 

Qsn+iis) = Q3„+4(l,s). 



APPENDIX D: CONSTRAINTS ON THE PADE 
COEFFICIENTS 

In this Appendix we consider the determination of the 
Pade coefficients {aj} and {bj} by application of the 
physical constraints (|5.7p or, equivalently, (|4.1ip . Sub- 
stitution of Eq. (|5.4p into Eq. (|4.7p . with the term e'' 
expanded in power series, yields 



M 



oo /j-2 



J2a,s^+' =Y^Qs^+Y^ Y.l,-2-iCi s^, (Dl) 



J=0 



3=2 \l=0 



where Cj is the jth coefficient in the series expansion of 
SAr(s)e* - X7]Am{s), i.e., 



-ATyOj, 0<j<M, (D2) 



min(j,Af) 



{l~k) 



., j>M + l, (D3) 



and Gj and 7-,- are defined in Eq. (|5.2ip . From the analysis 
of Eq. (jPip one obtains 



Cj=0, 0<j<d~l, 



Cj = aj-d, d<j<d+l, 



(D4) 



(D5) 



J=0 



0, m > + 1. 



In particular, Eqs. (jD4p and (jDSj), together with the 
definitions (jD2p and (jD3p . allow one to express the co- 
efficients {bj} in terms of {flj} for the cases N = d and 

N = d + l, 



bj — Sj^d + At? 



inin(j, A/) 



ai, 



0<j<d, (D7) 



M 



= ai - 1 + A?7^ — —ai, 



1=0 



(D8) 



where in Eq. (|D7p we have taken into account that ao = 
1. Notice that if iV = d, then bd+i = and Eq. (IPS]) 
provides an expression for one unknown aj in terms of the 
other ones. Expression (jPZp can be verified by induction 
using (|D4p and the following algebraic identity. 



(-ip 



m — 1 

y 

^-—^ j\{m — i)\ to! 



(D9) 



APPENDIX E: SOLUTIONS OF THE STANDARD 
RFA FOR d = 1, 3, AND 5 



The explicit solutions of Eq. (|5.26p are listed here for 
hard d-sphere systems in d = 1, 3, and 5 within the 
standard RFA approach presented in Sec. |Vl 



Hard rods 



The one-dimensional case {d = 1) corresponds to n = 
— 1. The recursive relation (IA9P yields O^iit) — 1/t, so 
that, according to Eq. (I4.ip . 



(El) 



G{s) ^ / drgir)e-'''. 
Jo 



Therefore, sG{s) is the Laplace transform of g{r). Taking 
n = -1 in Eq. (|5T5)) with A = 1 [cf. Eq. ([CT]) ]. one 
simply gets 



sG{s) 



1 



[77 + (1 — ri)s]e'' — 77 



(E2) 



The Laplace transform sG{s) can be easily inverted to ob- 
tain the well-known radial distribution function for hard 
rods [ill. Thus, the standard RFA becomes exact for 
d=l. 



Taking into account that, as a consequence of Eq. (|4.10p . 
the first n factors jj with j = odd vanish, it is obvious 
that jnH-dDl]) provide a total of (3dH- l)/2 = 3n + 5 
equations without unknown jj factors. Therefore, in 
agreement with Eq. (|5.7p . this is the number of con- 
straints on the Fade coefficients {bj} and {aj}. 



2. Hard spheres 

In the case d = 3 or n = 0, Eq. (|4.ip becomes 



G(s) = / drrg(r)e-"^ 
Jo 



(E3) 



17 



and so G{s) represents the Laplace transform of rg{r). 
According to Eq. (|5.18|) . 



G{s) 



1 + ais 



s2 l-12r7[03(s) + ai02(s)]' 



(E4) 



where we have taken into account that A = —12 and 
oo = 1. Equation (j5.26p reduces to Di = ai — ai, whose 
solution is 



ai = 



1 + ^/2 

1 + 27? ' 

From here one can easily get 



(1 + 277)2 



(E5) 



(E6) 



Equations (|E3 p - (|E6P define Wertheim-Thiele's exact so- 
lution of the PY integral equation for hard spheres 



Equation ()5.26|) yields Di — ai and D3 + ''JqDi = 0, 
where, according to Eq. ()5.25|) . 70 = 02 — 132 • The phys- 
ical solutions are 



ai 



3 - 37/(1 + 1002) 
1 — 677 



(E9) 



1 + 227? + 787?^ + 24773 + (1 - 77)(677 - 1)^ 



where ^ = a/I + 18?? + 677^. The associated contact value 
and isothermal susceptibility are 



3. Hyperspheres in five dimensions 

For a fluid of hyperspheres in 0? = 5 (ri = 1) the defi- 
nition of G(s) is 



dr r(l + rs)g{r)e 



(E7) 



Since uq — 3 and A = 240, Eq. (|5.18[) becomes in this 
case 



Gis) 



3 + fliS + 025^ 



1 + 24077[3(/)5(s) + ai04(s) + a203(s)] ' 



(E8) 



5(1+) 



- 1 + 337? + 87772 + 6773 
607?(1 - 77)3 



(l-7?)2 



(Ell) 



[5(1 + 677 + 3772) -2(2 + 3r?)e] . (E12) 



It can be checked that Eq. (|E8[) yields the same structure 
factor as that given by the solution of the PY integral 
equation for d = 5 obtained by a different method ^ iHi] . 
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